#### Perform the main specification under different EDD thresholds

## Setup the environment

setwd("~/research/coffee-dataverse/src")

do.origspec <- T
source("load.R", encoding="iso-8859-1")

library(lfe)
library(splines)

## Perform regression for each threshold

eddthresholds <- c(28, 30, 32, 34)
models <- list()
for (ee in 1:4) {
    print(eddthresholds[ee])
    thresholds <- c(10, eddthresholds[ee])

    for (col in paste0('above', thresholds)) {
        df[, col] <- 0
        for (monthsuff in monthsuffs)
            df[, col] <- df[, col] + df[, paste0(col, monthsuff)]
    }
    df$gdd1000 <- (df[, paste0('above', thresholds[1])] - df[, paste0('above', thresholds[2])]) / 1000
    df$edd1000 <- df[, paste0('above', thresholds[2])] / 1000

    mod <- felm(logyield ~ tmin + gdd1000 + edd1000 + prcp + prcp2 + state:ns(year, df=3) | Munic�pio + year | 0 | Munic�pio + year, data=df)
    models[[ee]] <- mod
}

## Produce table

library(stargazer)

stargazer(models, add.lines=list(c('Municipality FE', rep('Yes', 4)),
                                 c('Year FE', rep('Yes', 4)),
                                 c('State Trends', rep('Cubic', 4))),
          dep.var.labels="Log Yields", column.labels=paste(eddthresholds, " C"),
          covariate.labels=weatherlabels,
          omit=':', df=F, float=F)
